function amplification_factor(k)


p = params;

% compute A(t) numerically
lin_soln = stab(k, p);

t = lin_soln.x;
A = max(abs(lin_soln.y));
dA = diff(log(A)) ./ diff(t);

% [c, dcdt] = deval(lin_soln, lin_soln.x, p.N);
% dA = dcdt ./ c;


% computing s(t)
b_soln = base(p);
t = b_soln.x;
h0 = b_soln.y(p.N+1,:);

k_star = k * h0;
Ma_star = (h0 - 1 + p.beta) / p.beta * p.Ma;

for i = 1:length(k_star);
    
    p.Ma = Ma_star(i);
    s(i) = comp_eigs(k_star(i), p, h0(i)) ./ h0(i)^2;
    
end


subplot(3,1,1);
plot(lin_soln.x(1:end-1), dA, 'k', t, s, 'r--', 'linewidth',1);
xlabel('$t$','interpreter','latex','fontsize',12);
ylabel('Growth rate','interpreter','latex','fontsize',12);
% l = legend('$d(\log\, A(t)) / dt$','$s(t)$','location','northeast');
% set(l,'interpreter','latex','fontsize',11);
% legend boxoff
% ylim([-0.1 0.1]);

subplot(3,1,2);
dA = interp1(lin_soln.x(1:end-1), dA, t);
plot(t, abs(dA - s) / p.delta, 'k', 'linewidth',1);
xlabel('$t$','interpreter','latex','fontsize',12);
ylabel('Relative error','interpreter','latex','fontsize',12);
l = legend('$|d(\log A) / dt - s| / \delta$','location','southeast');
set(l,'interpreter','latex','fontsize',11);
legend boxoff

subplot(3,1,3);

S = exp(cumtrapz(t, s));
plot(lin_soln.x, A, 'k', t, S, 'k--', 'linewidth',1);
xlabel('$t$','interpreter','latex','fontsize',12);
ylabel('Amplification', 'interpreter','latex','fontsize',12);
l = legend('$A(t)$','$S(t)$','location','northeast');
set(l,'interpreter','latex','fontsize',11);
legend boxoff